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^ Abstract. We perform a computationl study of front speeds of G-equation models in time depen- 

1— I dent cellular flows. The G-equations arise in premixed turbulent combustion, and are Hamilton- 

^ Jacobi type level set partial differential equations (PDEs). The curvature-strain G equations are 

• also non-convex with degenerate diffusion. The computation is based on monotone finite differ- 

.^H ence discretization and weighted essentially nonoscillatory (WENO) methods. We found that the 
large time front speeds lock into the frequency of time periodic cellular flows in curvature- strain 
G-equations similar to what occurs in the basic inviscid G-equation. However, such frequency 

^ locking phenomenon disappears in viscous G-equation, and in the inviscid G-equation if time pe- 

^ riodic oscillation of the cellular flow is replaced by time stochastic oscillation. 
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>< 1. Introduction 

Front propagation in turbulent combustion is a nonlinear multiscale dynamical process [11, 14, 15]. 
In gasoline engine, the burning velocity determines the engine efficiency and the combustion cy- 
cle. Hence to measure and study the flame propagation speed is a fundamental issue in combustion 
experiments and theory. The comprehensive governing equations for flame propagation involve 
Navier-Stokes equations coupled with transport equations, which express laws of the fluid dynam- 
ics and the chemical reactions respectively. Simplified geometric models, such as the level set 
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G-equations [13, 11], are often more efficient in improving our understanding of such complex 
phenomenon. 

In level set approach [10], the so-called G-equation in turbulent combustion is: 

Gt + V(:sc,t)-DG + SL\DG\^0, (1.1) 

which describes the motion of the flame front at a prescribed flow velocity V(x, t) and a constant 
speed sl along the normal direction of the flame front called the "laminar flame speed". When the 
flame front is planar and the flow velocity is at rest, the front propagates at the speed sl- If the 
flow is in motion, the front is wrinkled in time and eventually propagates at an asymptotic speed 
St in each specified direction, which is called the "turbulent flame speed". Through the framework 
of homogenization theory, the existence of the turbulent flame speed of (1.1) has been rigorously 
established for periodic flows [16, 3] and random flows [7, 4]. 

As fluid turbulence is known to cause corrugations of flames, additional modeling terms are of- 
ten incorporated into the basic G-equation (1. 1) on physical grounds [11]. We shall study turbulent 
flame speeds characterized by such extended G-equation models involving curvature, diffusion or 
strain effects. The curvature G-equation is: 

Gt + V{^,t)-DG + sl\DG\ = d\DG\diY (^-^^ , (1.2) 

which comes from adding mean curvature term to the basic motion law. If the curvature terms is 
Unearized, we arrive at the viscous G-equation: 

Gt + V{^,t)-DG + SL\DG\=dAG, (1.3) 

which is also a model for understanding the numerical diffusion [10]. The strain G-equation is: 

Gt + y{^,t)-DG+(^SL+d^^j^^^^ \DG\ = d\DG\diy (J^^ . (1.4) 

The strain term n-DVn comes from the flame surface stretching rate. We will give a brief derivation 
of all G-equation models later. 

Our goal is to study turbulent flame speed st and its dependence on the velocity field l^(x, t). 
Cellular flows are two dimensional vortical flows with nontrivial geometric streamline structures, 
hence often adopted for mathematical and computational study of st- Let us first consider the 
steady cellular flow: 

V(x, t) — A-{— sin(27rx) cos(27r7/), cos(27ra;) sin(27r7/)). 

which is a two dimensional periodic incompressible flow with parameter A as the flow intensity. 
The streamline of the cellular flow consist of a periodic array of hyperbolic (saddle points, sepa- 
ratrices) and elliptic (vortical) regions. For inviscid G-equation, [9, 1] showed that st is enhanced 
by the cellular flow with growth rate st = 0{A/ log A). For viscous G-equation, it is proved in 
[5] that the diffusion term causes strong front speed bending (saturation), and the growth of st 
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slows down to uniform bound: st — 0{1). For curvature and strain G-equations, although the an- 
alytical justifications are still lacking, numerical results [6] showed that the curvature term causes 
weaker front speed bending than the regular diffusion term. Also the strain rate term causes flame 
quenching, that is, st starts to decrease for larger A and eventually drops to zero. 

Next, we consider the unsteady (time dependent) cellular flow. For inviscid G-equation (1.1) 
and time periodic cellular flow, it is observed in [2] that the front propagation may be synchronized 
by the spatial and temporal periodicity, and st is approximately frequency times a rational number. 
Such phenomenon is called "frequency locking". The time periodic and random cellular flows are 
also considered in [8] within a reaction-diffusion-advection (RDA) equation model where tempo- 
rally stochastic perturbations of the cellular flow also cause strong resonance (although the number 
of resonance peaks reduces to one). 

In this paper, we study turbulent flame speeds in unsteady cellular flows for extended G- 
equation models (1.2), (1.3), (1.4). The organization is as follows. In section 2, we give an overview 
of G-equation models and refer the associated numerical schemes based on monotone finite dif- 
ference methods to [6]. In section 3, we show numerical findings. For time periodically shifted 
cellular flows, we found that frequency locking of the front propagation persists in curvature and 
strain G-equations, and the locking phenomenon is robust with respect to the coefficient d of the 
curvature or strain term. On the other hand, it disappears in the viscous G-equation (1.3), and 
in time random cellular flows where st decreases in the oscillation "frequency" (the reciprocal of 
correlation length). For time periodically amplified cellular flows, all results are similar except that 
if the amplitude oscillation frequency is too high then the front becomes insensitive to the cellular 
flow. We end with concluding remarks in section 4, and acknowledgements in section 5. 



2. G-equation Models 

In the thin reaction zone regime and the corrugated flamelet regime of premixed turbulent combus- 
tion [11], the flame front is modeled by a level set function {G'(x, t) = 0}, which is the interface 
between the unbumed fuel {G < 0} and burned fuel {G > 0}. The trajectory of a particle x{t) on 
the interface is given by a velocity field and a laminar speed: 

dx 

— = V(x,t) + Sin, (2.1) 

where the positive constant sz, is called the laminar flame speed and n = DG/\DG\ is the normal 
velocity (D: spatial gradient). In terms of G, the motion law (2.1) gives the inviscid G-equation 
(1.1). To take account of the effect of flame stretching, the surface stretch rate may be added as a 
first order correction term on laminar flame speed: 

sl = sl - d{K + S), 

K = div(n), S = —n-DV-n. 

Here d is called the Markstein diffusive number, k is the mean curvature and S is called the strain 
rate. Replacing s^, by sl in the motion law (2.1), we obtain the strain G-equation (1.4). If we 
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include the motion by curvature only, we have the curvature G-equation (1.2). If the curvature 
term is further linearized, we arrive at the viscous G-equation (1.3). 

Suppose the flame front propagates in x-direction in two dimensional space (x = {x,y)). 
Consider the stripe domain IR x [0, 1], and the burned region at time t is {x G M x [0, 1] : G(x, t) < 
0}. Denote A{t) the area that the burned region has invaded during time interval (0, t): 



(X' indicator function.) We refer A'{t) as the instantaneous burning velocity. Then turbulent flame 
speed is the large time Unear growth rate of A{t) or the average of A'{t): 



In numerical simulation, consider G-equation with initial data G(x, 0) = x. If V(x, t) is spa- 
tially periodic, we can reduce the computational domain to [0, 1]^ by imposing the affine periodic 
condition: 



Based on the framework of finite difference computation of Hamilton- Jacobi equations [10], the 
first order derivatives (flow velocity, laminar speed and strain terms) are discretized as mono- 
tone numerical Hamiltonian and approximated by weighted essentially nonoscillatory (WENO) 
scheme. The second order derivatives (curvature and diffusion terms) are approximated by cen- 
tral differencing. For the explicit time step discretization, total variation diminishing Runge-Kutta 
(TVD RK) scheme is used, and the time step size restriction is given by the CFL condition. See 
[6] for details of the numerical schemes. We compute (2.3) with spatial grid size 200 x 200. 



3. Numerical Results 

Our first set of numerical results is on the time periodically shifted cellular flow: 

V(x, t) = A-{— sin(27ra; + B sm{2nut)) cos{2ny), cos(27ra; + B sm{2nujt)) sm(2Tcy)). 

The flow oscillates temporally in the x-direction with amplitude B and frequency cu. We fix = 1, 
A = A, B = 1 and choose different values for d to see the effect of the curvature, diffusion and 
strain terms. Note that if d — then all G-equations are identical to inviscid G-equation, and we 
choose d small enough to avoid the strain term overpowering sl- 

Figure 1 shows the plot of A'{t) for inviscid G-equation. For uj — 2.8, A'{t) evolves into a 
periodic function after a transient time. In this case frequency locking occurs, and we can simply 
compute St by taking the average of A'{t) over a few periodic intervals. For ui = 3.2, A'{t) 





(2.2) 




(2.3) 
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Figure 1 : Plots of A'{t) for inviscid G-equation with time periodically shifted cellular flow, s l 
A = A,B = 1. Right: u = 2.8. Left: u = 3.2. 



1, 



remains chaotic over time, and st is evaluated by taking average of A'{t) over a large time interval 
as (2.2). 

The frequency locking effect can be understood as follows. In the transient time, the front 
propagation structure is synchronized with the temporal and spatial periodicity of the cellular flow 
and evolves into a periodic steady state. During a multiple of the oscillatory period At = M ■ l/cu, 
the flame front translates through a multiple of the cell AA{t) = N. {N, M: integers.) Then the 
turbulent flame speed is: 

st= . = tt-^- (3-1) 

At M 

Figure 2 shows the plot of turbulent flame speed as a function of the oscillatory frequency st = 
sxiuj) for all G-equation models with various d. For inviscid, curvature and strain G-equations 
(1.1),(1.2),(1.4), frequency locking is observed. Here sr(cu) is piecewise linear with rational slopes 
r = st/uj- In each linear segment the front propagation locks into the same periodicity pattern, 
that is, the same A^, M in (3.1) and hence same r. Between the linear segments, the locking effect 
loses stability and st{co) decreases. Notice that frequency locking phenomenon is also robust with 
respect to d, that is, st may be the same even as d varies (different curves may coalesce on bold 
line segments of top and middle panels of Figure 2). This is different from that in steady cellular 
flow where st is strictly decreasing as d increases [6]. 

For viscous G-equation, however, no frequency locking appears. Instead the dissipation effect 
causes the level set function to converge to the stationary solution: 

(j(x, t) = —Sxt + X + ^(x, t) 

with u{x,t) a periodic function in [0, 1]^ x [0, l/co]. In the plot of st{uj), we see two resonant 
peaks for smaller d, then the peaks start to merge and disappear for larger d. 

Our second set of numerical results is when V is the time periodically amplified cellular flow: 

V(x, t) = y4 sin(27rci;t) ■ (— sin(27rx) cos(27r?/), cos(27rx) sin(27ry)). 

Figure 3 is the plot st = st{uj) for all G-equation models and various d. When cu is small, all 
results are similar to the time periodically shifted cellular flow. When cu is larger, however, the 



5 



Y. Liu et al. 



Hirbulent Flame Speeds of G-equation Models in Unsteady Cellular Flows 




Figure 2: Turbulent flame speeds of G-equation models with time periodically shifted cellular 
flow. Above: inviscid/curvature G-equations. Middle: inviscid/strain G-equations. Below: invis- 
cid/viscous G-equations. The bold lines indicate the frequency locking pattern r — st/oj. 
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direction of the flow changes so fast in time that the flame front cannot be really wrinkled, and 
so remains almost planar. Hence st{uj) is not much responsive to the cellular flow patterns, and 

decreases to s/, as — )■ +00 in all the G-equations. 

Finally we consider the time stochastically shifted cellular flow: 

V(x, t) — A-{— sin(27rx + ^{t)) cos(27ry), cos(27rx + ^{t}) sin(27ry)). 

Here 7(t) is the Omstein-Uhlenbeck process defined by the Ito equation: 

d-f{t) = -a-f{t)dt + /3dW{t), 

where W{t) is the Wiener process, a and /3 are positive constants. We choose /3 — \plo?l'^ so 
that the power spectrum is invariant in «, the reciprocal of the correlation length of 7(i). Figure 
4 shows the plot s^ip) for inviscid G-equation. We see that sj- always decreases as ol increases, 
stochasticity in time kills the resonance phenomenon. This is different to that in RDA equation 
model where one wider resonance peak still remains under time stochastic perturbation of cellular 
flows [8]. 



4. Concluding Remarks 

We have studied numerically the basic and the curvature-strain dependent G-equation models and 
their corresponding turbulent flame speeds in various unsteady cellular flows. Numerical results in- 
dicate that the frequency locking phenomenon occurs in inviscid, curvature and strain G-equations 
but not in viscous G-equations. It also disappears in time stochastically shifted cellular flows. In 
the future, we plan to study more general unsteady (time periodic or random) flows, including three 
dimensional flows for G-equation models. 
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Figure 4: Turbulent flame speeds of inviscid G-equation with Omstein-Uhlenbeck process being 
time oscillatory shift in cellular flow, as a function of a, the reciprocal of the correlation length, or 
"oscillation frequency". 
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